Pseudochaos and low-frequency percolation scaling for turbulent diffusion in 

magnetized plasma 
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The basic physics properties and simplified model descriptions of the paradigmatic "percolation" 
transport in low-frequency, electrostatic (anisotropic magnetic) turbulence are theoretically ana- 
lyzed. The key problem being addressed is the scaling of the turbulent diffusion coefficient with 
the fluctuation strength in the limit of slow fluctuation frequencies (large Kubo numbers). In this 
limit, the transport is found to exhibit pseudochaotic, rather than simply chaotic, properties asso- 
ciated with the vanishing Kolmogorov- Sinai entropy and anomalously slow mixing of phase space 
trajectories. Based on a simple random walk model, we find the low-frequency, percolation scaling 
of the turbulent diffusion coefficient to be given by D/uj oc Q 2 ^ 3 (here Q ^> 1 is the Kubo number 
and u) is the characteristic fluctuation frequency). When the pseudochaotic property is relaxed the 
percolation scaling is shown to cross over to Bohm scaling. The features of turbulent transport in 
the pseudochaotic regime are described statistically in terms of a time fractional diffusion equation 
with the fractional derivative in the Caputo sense. Additional physics effects associated with finite 
particle inertia are considered. 



I. INTRODUCTION 

It is generally agreed that the presence of low- 
frequency, long-wavelength fluctuations in hot magne- 
tized plasma may have a deteriorating effect on the 
plasma confinement leading to an anomalously high heat 
and energy transfer across magnetic field lines as com- 
pared to purely collisional values, a phenomenon known 
as anomalous or turbulent transport. While a first prin- 
ciples theory of turbulent transport is not at hand, the 
practical evaluation of the transport level is often ob- 
tained by scaling relations. An important problem inves- 
tigates the scaling laws for the diffusion coefficient as a 
function of the fluctuation strength or the so-called Kubo 
number Q ~ u± / lu£;±_ , where uj_ is the characteristic flow 
velocity, £x is the cross-field correlation length, and lo is 
the typical fluctuation frequency. 

Since the early studies of Isichenko and co-workers 
1,13, Jj, Jj| it has been discussed by a few authors 
5, 111, 0, S, S E3] that the diffusion coefficient due to 
turbulence exhibits a power-law dependence D/oj oc Q 7 , 
where the scaling is quasilinear-like (7 = 2) for Q -C 1 
and percolation- like (7 = 7/10) for Q 3> 1 (as distinct 
from Bohm scaling with 7=1 11]). It has also been 
discussed [f| [IH, [l3[ that the percolation scaling with 
7 = 7/10 is however not exact and that the true value of 
7 is actually smaller than (although remarkably close to) 
Isichenko's original estimate, 7/10. More so, an improved 
percolation scaling with 7 = 2/3 12] has been proposed 
in connection with a fractional generalization [3, 0, [l6| 
of the Fokker-Planck-Kolmogorov equation. This devia- 
tion from the earlier predicted 7 = 7/10 finds support 
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in the numerical simulation of anisotropic magnetic field 
turbulence [9I [l3j. 

Although small, the observed discrepancy needs to be 
addressed. Indeed this discrepancy leads to noticeable 
variation in the predicted transport level (because of the 
large Q ^> 1). Apart from the numerical differences, the 
basic physics origin of the percolation scaling, as well as 
the fundamental reason for the deviation from 7/10, has 
not been clearly understood. 

In this work, we expose a few crucial physics aspects 
behind the percolation scaling of the turbulent diffusion 
coefficient. We find that a mathematically consistent ap- 
proach to turbulent diffusion in the limit of low frequen- 
cies (large Kubo numbers) can be obtained within the 
concept of pseudochaos (random non-chaotic dynamics 
with zero Lyapunov exponents) [l7], [H, [H| . Our analysis 
displays a few characteristic features of pseudochaos dif- 
ferentiating it from the more intuitive, chaotic behavior. 
We confirm the percolation scaling with 7 = 2/3, basing 
our considerations on a simple random-walk model in 
fractal geometry. Yet, a slightly smaller value 7 ~ 0.66 
could be advocated involving subtleties of the random 
walks at percolation. 

More so, we demonstrate that the diffusion on per- 
colation systems is described by a non-Markovian dif- 
fusion equation, with the non-Markovian property ex- 
pressible in terms of a fractional time derivative in the 
Caputo sense [20| . Here the stress is on "Caputo" and 
on subtleties of definition of the fractional derivative 
operator. Being conceptually very similar, this defini- 
tion is somewhat different from the definition in the 
Riemann-Liouville sense [2l| . of wider use in the appli- 
cations. Indeed fractional diffusion equations with the 
Riemann-Liouville derivative have, in the various con- 
texts, been studied and discussed in the literature (Refs. 
[H El H| for reviews) . 

The paper is organized as follows. The concepts of 
chaos, pseudochaos, and resonances in Hamiltonian dy- 
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namics are discussed first (Sec. II), followed by a deriva- 
tion of the percolation scaling from a random walk model 
in fractal geometry (Sec. III). Fractional derivative equa- 
tions are obtained in Sec. IV in the framework of gener- 
alized memory function formalism. Next issues related to 
particle inertia are considered (Sec. V). We summarize 
our findings in Sec. VI. Applications of this study pertain 
to both geo-space [26|, [27j and fusion [H, [2^| plasma, and 
by mathematical analogy to problems outside the plasma 
physics. 

II. HAMILTONIAN, RESONANCES, AND 
PSEUDOCHAOS 

To begin, we formulate the transport problem as the 
Hamiltonian problem 

dx _ dH(x,y,t) dy _ dH(x,y,t) 

dt ~ dy ' dt ~ dx U 

where x and y are coordinates in the plane perpendicular 
to the magnetic field, and H(x,y,t) is a time-depending 
Hamiltonian. 



A. Hamiltonian problem - examples 

The physics included in Eqs. (p} encompasses, beside 
others, the following model realizations: 

(i) Field-line diffusion by anisotropic magnetic turbu- 
lence. Early work on this problem is due to Rosenbluth 
et al. 30]. The magnetic field model is set as a sum of 
a constant homogeneous background field Bo = Bqz in a 
straight-cylinder geometry (here z is a unit vector in the 
axial direction) and a static transverse magnetic pertur- 
bation, <5Bj_(r_L, z), small compared to the background 
field: 8B±/Bq <C 1. The equations for the magnetic 
field lines are cast in the Hamiltonian form by defin- 
ing 5ti± = VAu(r±,z) x z, where rj_ = (x,y) is the 
position vector in the perpendicular plane, • z = 0, 
V = d/dr±, H(r±,z) — A\i{t±,z)/Bq is the Hamilto- 
nian, and the axial coordinate z is considered as "time." 
The flow velocity is given by uj_ = VH(r±,z) x z = 
S~B±/Bq 7 making it possible to evaluate Q ~ u±/lu^j_ ~ 
(l/2ir)(6B±/Bo)(£u/£±), where SB^ is the rms magnetic 
fluctuation. In writing the last form we took into account 
that the "frequency" of the field variation due to the z 
dependence is u> ~ 27r/£ii, where £|| is the typical par- 
allel correlation length. It is noted that Q > 1 implies 

3> 1 with a large margin, and hence, it is required 
the magnetic turbulence be strongly anisotropic. 

(ii) Diffusion of guiding centers by the Ej_ x B drift. 
This problem refers to the well-known problem of particle 
diffusion by low-frequency, electrostatic, microturbulcnce 
[S HH H3- This is usually one of the most important 
transport problems in a fusion plasma. The Hamiltonian 
is written as H(r±,t) = —(l/B )$(r±,t), where <&(rj_, t) 
is the electrostatic potential in the perpendicular plane, 



such that Ei(ri,t) = — V$(rj_, i), and Bq is the ampli- 
tude of a confining, uniform magnetic field. The control 
parameter is defined by Q — Ej_/lu£±Bq, where E± is 
the rms electric field, and u is the typical fluctuation 
frequency. 

(iii) Diffusion of guiding centers by the VB drift. 
Consider the magnetic field in the z direction and sup- 
pose that it is spatially nonuniform in the perpendicu- 
lar plane in accordance with B = B(r±,t)z, B(r±,t) = 
B + 5B(r±,t), and SB/B <C 1. The time dependence 
in 5B(r±,t) is assumed to be very slow, so that the cor- 
responding variation frequency is, with a large margin, 
small compared to a characteristic cyclotron frequency, 
i.e., LojiQ c -C 1. It is understood that the time vary- 
ing magnetic perturbation generates a time and spatially 
varying electric field because of Faraday's law. A con- 
sistent description of the guiding center motion, then, 
will require that the Ex X B and polarization drifts are 
taken into consideration. A simplified yet relevant model 
is obtained for situations in which the spatial inhomo- 
geneities of the magnetic field are attributed a major 
role in determining the characteristic transport proper- 
ties of the considered system as compared to the explicit 
time dependencies of the magnetic fluctuation and induc- 
tive electric fields. For instance, such model conditions 
can occur for the transport of solar wind-like plasma in- 
side the magnetosphere as a consequence of the internal 
fine structure of the Earth's magnetopause [2(|. In this 
spirit, if we assume that the VB drift velocity is the 
dominating velocity in the Alfven approximation, and 
neglect the higher-order corrections due to the inductive 
electric field, we can write u_l — TMb(z x VB)/qB . 
where the upper sign corresponds to a negative charge 
and we have, following Refs. [H, l34j . introduced the 
first adiabatic invariant [Ib — mv\/2B (here v±_ is the 
perpendicular particle velocity, m is the particle mass, 
and other notations are standard). In order for the 
drift approximation be valid, we also need to require 
(vj_/lj c )\VB/Bo\ <C 1. We note in passing that the VB 
drift is a finite Larmor radius effect. The Hamiltonian is 
defined by H(r±,t) = ±{n B /qB )SB(r ± ,t). This model 
Hamiltonian isolates the effect of the magnetic inhomo- 
geneity. It includes the magnetic gradients slowly chang- 
ing with time. Effects due to the time varying electric 
field, left apart in the present model, will be considered 
in more detail in Sec. V. 



B. Equipotentials - fractal structure 

At a given time t the equipotentials of the Hamiltonian 
are defined by y — y(x, H,t — const). By standard rules 
the area embraced by a closed equipotential is obtained as 
I = § ydx (this being the action of the Hamiltonian sys- 
tem), while the frequency of the orbital motion is given 
by Q(I) = dH{I)/dI ~ \VH(r ±: t)\/\VI(r ± ,t)\ ~ u±/t, 
where I ~ |V/(rj_, t)\ is the circumference of the equipo- 
tential contour, and t is kept as a fixed parameter. As- 
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suming circular equipotentials one writes I oc r\ , £ oc r± , 
and I oc £ 2 . In the case of strongly shaped equipotential 
contours, however, the latter scaling relations may not 
hold true. Of particular interest are the equipotentials 
characterized by self-similar geometry in a broad range of 
spatial scales best described as fractals [35[ . The implica- 
tion of such equipotentials in diffusion by low-frequency 
turbulence is addressed shortly. For fractal equipoten- 
tial contours the simple scaling I oc £ 2 is generalized to 
j x g2/d h (-^g so _ ca j] ec ] "area-perimeter relation" [361 ]) 
where dh is the fractal dimension of the equipotential 
line and we have suppressed the normalization param- 
eter for simplicity. The fractal dimension dh is defined 
by the scaling £ oc r± (instead of £ oc rj_ in non-fractal 
geometry). The dh values are generally fractional larger 
than 1. Eliminating £ by means of the area-perimeter 
relation one obtains £1(1) oc uj_/I d '^ 2 . 



scales shorter than a and that there is a transition to a 
statistically homogeneous distribution of the fluctuations 
at length scales longer than £x- 

Consistent with the above assumptions regarding the 
guiding center approximation, we shall also assume that 
a is large compared with the typical scales of the parti- 
cle gyro motion. When the particle gyro radius is not 
small enough to match with this condition, then the par- 
ticle motion can be considerably modified by the multi- 
scale geometry of the turbulent field. In some cases the 
particle trajectory can be extremely complicated. One 
example of this complex motion is the case of "strange 
diffusion," discussed in Refs. [HI, H3]. We also note that, 
when £j_/ a — > 00 > the expansion in Eq. @ reproduces 
the so-called fractional Brownian function, which offers a 
particularly clear example of fractal behavior [3^, H2, 53] . 



C. Hamiltonian - Fourier representation 

In the discussion so far we have not made any specific 
assumptions about the structure of the H field in space 
and time. We now need to be more detailed in connec- 
tion with the notion of "turbulence." It is convenient to 
think of turbulence as of collection of plane waves, with 
amplitudes 6 nj k and random phases </> n ,k, making it pos- 
sible to expand the Hamiltonian into the Fourier series to 
give (with k being the wave vector in the perpendicular 
plane) 

H(r±,t) ~y^y^ b n<k cos(k • r ± - uj n t + n ,k). (2) 

n k 

Next, it is generally agreed that there is a broad, isotropic 
wave vector spectrum characterized by a power-law. In 
2D geometry we can write 

P(k)^2KkY,b 2 nM ^k- a . (3) 

n 

Actually a more precise way of looking at the spectral 
characteristics of the turbulence would be to introduce a 
fully 3D spectrum, squeezed to a "pancake" in the per- 
pendicular plane. Such spectra are considered in Ref. 

The next point to be addressed is the a value in Eq. ((3| . 
For the present analysis, we take 1 < a < 3. Together 
with the inverse power-law energy density distribution, 
the latter condition implies that the fluctuations are self- 
similar, in the sense of Refs. [3511361]. Indeed this prop- 
erty of self-similaraty finds suppor t in the many different 
experiments [13, H3, 1H, [H, So, |4l|] . From the standpoint 
of formal treatment, the fluctuations being self-similar 
pave the way for the application of fractal geometry as 
proposed in the model below. 

To match with the physics natural limitations, the 
spectrum in Eq. ([3]) is characterized by finite cutoffs on 
both sides, i.e., a <C 2ir/k <C Without loss in gener- 
ality, we shall assume that all fields are smooth at length 



D. Percolation property and resonances 

Because of the sign symmetry of the Hamiltonian in 
Eq. the zero-set H(r±,t = const) = contains a per- 
colating equipotential line [J, HH , which, by its geomet- 
ric meaning, stretches the entire system. The percolat- 
ing line is the channel through which turbulent transport 
penetrates to the large scales [1] . One way to describe the 
transport is to directly link it to the complex folding of 
the percolating line in the real space and in this fashion to 
predict the anomalous scaling laws for the turbulent dif- 
fusion coefficient as suggested in the early work [l| . This 
approach, however, leads to some difficulties associated 
with the fractal dimension of the transport process as it 
basically ignores the dynamical features of the transport 
[l2l ]. A correct way to deal with the properties of tur- 
bulent diffusion is to account for resonances between the 
orbital motions of the particles and the time variation of 
the Hamiltonian. As is well known, when phase space 
trajectories are subject to local instability, the matching 
fluctuation and orbital frequencies or the matching har- 
monics of these favor departure from the exact periodic, 
energy conserving motion thus giving rise to transport 
phenomena in phase space (Refs. (4j| [4(| [47| for the 
full discussion). In the limit of very low frequencies, the 
conditions for a resonance require the corresponding ex- 
cursion periods to diverge. Mathematically, this can be 
satisfied for fractal equipotential contours, characterized 
by the diverging lengths due to the structure on many 
scales. Analysis of Ref. [44| have shown that the equipo- 
tentials of the H field in Eq. @ occurring in close prox- 
imity to the percolating line are indeed described by the 
fractal structure in the limit £j_ /a — > 00 and their fractal 
dimension was found to be dh — 1.32 ±0.01. This fractal 
dimension is close to the value dh = 1.37±0.03, obtained 
in Ref. [48|, and is clearly smaller than the hull exponent, 
dh = 1.75, used by Isichenko and co-workers [3, @, 0]. 
Anticipating our result in Sec. Ill, here we mention that 
the "percolation" scaling of the turbulent diffusion coef- 
ficient won't however involve the dh parameter directly, 
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in contrast with the conclusion of Refs. [H, @, 0). 



E. From the percolation property to pseudochaos 

When the fluctuations are very slow, we expect the 
resonance conditions to be naturally satisfied in vicin- 
ity of the percolating line, provided that the number of 
modes in Eq. ([2]) is large enough. As is already men- 
tioned in the above discussion, this line, being strongly 
shaped, is characterized by the fractal geometry in the 
limit £j_/a — > oo. Fractality, in its turn, guarantees the 
existence of resonant orbits for tu n — > 0. It is noted that, 
when uj n — > 0, the resonances strongly overlap in a very 
narrow layer enveloping the percolating line. This obser- 
vation actually refers to the general properties of Hamil- 
tonian dynamics near the separatrices, discussed in Ref. 
pTi ] . It is the overlap of resonances that introduces insta- 
bility into the motion. If one defines the characteristic 
frequency, u>, one finds the width of the resonance layer 
to be Slu cx ^/uj (as compared to the distance between the 
resonances Alu cx lo). Note that Slu ^S> Aw for lo — ► 0, 
and that the width of the resonance layer shrinks to zero 
(as y/uj), while the density of the resonances in the layer, 
defined as inverse distance Aw, diverges as 1/lu. Inside 
the layer, the dynamics are basically random because of 
the many resonances present [HI, H3| ■ On the other hand, 
the concentration of the resonance properties on a fractal 
structure, supported by £j_/a — > 00, strongly reduces the 
phase space available for the random motions. As a con- 
sequence, the chaotic trajectories separate anomalously 
slowly. "Anomalously slowly" means sub-exponential for 
lo — > (see the discussion below), which is equivalent to 
saying that the Lyapunov exponents, thought of as func- 
tions of position in phase space, vanish at percolation. 
In view of the vanishing Slu — > the phenomenon can be 
envisaged as "stickiness" to the percolating line. More 
so, the vanishing Lyapunov exponents, combined with 
the inherently random character of the motion, indicate 
that the dynamics are pseudochaotic, in the sense of Refs. 
[I?], [H, [l9|, instead of being just chaotic. This is the key 
point. Based on the notion of pseudochaos one al- 
most straightforwardly derives the "percolation" scaling 
of the turbulent diffusion coefficient in the limit of slow 
fluctuation frequencies. 



F. K parameter 

In order to make these intuitive arguments more pre- 
cise, we proceed as follows. In the basic theory of Hamil- 
tonian chaos one writes the number of overlapping reso- 
nances as y/K, where the K value (the so-called nonlin- 
earity parameter) is obtained from [45l |46l [47j 



Consistent with the above orderings, \[K cx ^ 1 in 

the resonance layer. Also, if we denote the Lyapunov ex- 
ponents of the Hamiltonian system in Eqs. (p} by ±cr, we 
can order \a\ ~ /iks — l/ r c- Here /iks is the Kolmogorov- 
Sinai entropy and t c is the characteristic mixing time. 
For K 3> 1, the mixing time is written as r c ~ 1/2lo\d.K 
(Refs. [45|,|47[ for the full discussion). When lo — ► 0, the 
K parameter goes to infinity as 1/lo, while the mixing 
time diverges as r c cx — I/2wlnw — > 00. Accordingly, 
|er| ~ /iks — — 2lo\tilo — > at percolation. Thus, while 
the dynamics arc intrinsically random due to the many 
overlapping resonances, the available phase space is in- 
sufficient to host the chaos. Indeed the mixing happens 
to be sub-exponential permitting the system exhibit only 
a pseudochaotic behavior. 

Let us now rewrite Eq. |4]) in a more insightful form. 
Substituting £1(1) cx uj_/I dh ^ 2 and performing the triv- 
ial differentiation over / leads to if cx u±/lo£ 7 where the 
area-perimeter relation I cx £ 2 / dh has been applied and 
we have omitted the numerical constant factor for sim- 
plicity. For the fractal equipotential contours near the 
percolating line, the circumference £ cx ^ , making it 
possible to evaluate K cx u±/lo^\ Remembering that 
the Kubo number Q cx it is found that 



Q/K cx d 



(5) 



K = — 

LO 



dft(I) 



dl 



(4) 



Hence, in the pseudochaotic regime, Q ^> K due to frac- 
tality (£±/a > 1, dh > 1). When d h — > 1, the Q and K 
parameters have the same order. 



G. Kolmogorov-Sinai entropy 

Eliminating K by means of Eq. ((5]) we can also express 
the Kolmogorov-Sinai entropy as a function of the Kubo 
number as /iks — 2ljItiQ — 2uj(dh — l)m£j_ — 2wlnC, 
where the last term in the sequence comes from the 
(suppressed) normalization constant factor in the above 
scaling relation. It is noted that, to the leading order, 
Hks is a logarithmic function of Q. More generally, this 
logarithmic behavior can be slowed down by dh nonlin- 
ear ly depending on the Kubo number. The implication is 
that fractal geometries of the contours on which the reso- 
nances occur can differ with the control parameter. This 
effect will be discussed elsewhere [4!| . Our result for /iks 
is in contrast with Isichenko's finding /iks — wQ 1//2 ln<5 
for Q ^> 1 (Eq. (4.55) in Ref. 0) and the square-root- 
like scaling versus Q, argued to be "universal." Indeed it 
is found in direct numerical computations [49 . [50] that 
the square-root-like scaling of the entropy is not repro- 
duced and that the behavior is actually slower than a 
power-law, tending to saturation. We interpret these re- 
sults as consistent with implications of the slow mixing 
and pseudochaos. 
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H. Summary 

The end result of the above reasoning is that turbu- 
lent diffusion in the limit of very low frequencies is char- 
acterized by the random dynamics squeezed to a sub- 
set of phase space with fractal geometry. In this regime 
Q K ^> 1, whereas for the classical chaotic behav- 
ior [i^, [47|, characterized by the wide domains of ran- 
dom motion, Q ~ K 3> 1. Overall, the dynamics bear 
signatures enabling to associate them with pseudochaos 
(random non-chaotic dynamics with zero Lyapunov ex- 
ponents) [l7l | . 



III. PERCOLATION SCALING FOR 
TURBULENT DIFFUSION 

We are now in position to obtain the scaling law for the 
diffusion coefficient in the parameter range of slow fluctu- 
ation frequencies (large Kubo numbers). The calculation 
is based on the general scaling form 

{rl(t))=eAt/T*)f(t/T*), (6) 

where r* is the characteristic diffusion time to the dis- 
tance and / is a scaling function, such that /(oo) = 
const. The form in Eq. ^ is similar to that considered by 
Gefen et al. |o"l| for anomalous diffusion on percolation 
clusters (in their model r* oc where 9 is the index 

of anomalous diffusion), and earlier by Straley [52| . 

Based on the above discussion, we assume the / func- 
tion in Eq. ([6]) has already reached its asymptotic value 
after the characteristic mixing time t ~ t c . That is, after 
this time, the dispersion in Eq. is linear. On shorter 
time scales, the dispersion deviates from linear because of 
the concentration of the transport process on a fractal ge- 
ometry. Since the dynamics are nevertheless random, we 
can envisage them as a random walk process on a fractal 
cluster at percolation. This process has been discussed 
as a simple model process for diffusion in disordered sys- 
tems ("the ant in the fractal labyrinth") [H^,|5^{. In those 
models a particle (random walker) is assumed to hop in 
random manner between nearest neighbor sites of a frac- 
tal lattice moving a step of length a each a/uj_ seconds. 
In more advanced models the condition that the steps oc- 
cur at fixed time intervals is relaxed [55|. Note that, in all 
these models, the fractal lattice is assumed to be station- 
ary. In turbulent diffusion, this assumption holds true for 
only short time lags t < r c <C l/u>. At these short times, 
we can require the / function in Eq. ^ correspond with 
the scaling function in the simple random walk model 
|5ll ]: f(s) oc s~^M 2 + e ). Here /i and v are the percolation 
indices, whose definitions are explained in major reviews 
0) HI, |57| ■ It is remarked that t c has the sense of the 
characteristic cross-over time scale where the / function 
changes from the initial-time power-law behavior to the 
asymptotic constant value. The diffusion coefficient is 
evaluated by noting that /(oo) ~ (tc/t*)"^ 2 "^, then 



making use of Eq. ((6|) to obtain 

D ~ (£i/2T»)(r c /T,)^ M2+e) . (7) 

Focusing on the random walks at percolation we some- 
how expect the fractal range be wide enough, ideally 
£_i_/a — » 00. On the other hand, the slowness of the fluc- 
tuations implies that Q ~ uj_/lu£j_ — > 00 for to — ► 0. For 
the purpose of formal orderings, we require £j_/a — * 00 
diverge faster than Q — ► 00. That is, given a finite fluc- 
tuation frequency, we set £j_/a ^> u±/lo^±. The impli- 
cation is that the fluctuations should not be too slow 
for the actual (finite) size of the system, otherwise the 
dynamics retain a deterministic character. If we de- 
fine t* ~ £j_ / au± , from the last inequality we also get 
wt, ^> 1 showing that the diffusion time to the distance 
must be large compared to the period of the field. 
Because of this, we expect fractal correlations to decay 
before the particles have crossed the entire system. Thus, 
in our model, the "anomalous" scaling r* oc £ 2 _ +e [5l| is 
invalidated in the limit £j_ — > 00 as a consequence of 
the time varying Hamiltonian. Instead, the conventional 
scaling law for diffusion, r* oc applies. Remembering 
that t c ~ l/2uhxK, we can cast the diffusion coefficient 
in Eq. ([7]) in the form 

D~2-'(\nK) 1 -'(a/U)'^e±Q\ (8) 
where the exponent 7 is defined as 

7 =l-/VK2 + 0). (9) 

This is the desired result. We have D/oj oc Q 7 , where 
7 depends on the percolation indices, but the expression 
is different from the one obtained by Isichcnko [4]. Note 
that D/wis not a simple scaling function of to (because of 
the K parameter varying with frequency) . The correction 
is expressible as a power of the logarithmic dependence 
oc (— In a;) 1-7 , where to — > 0. 

A. General expression 

It is convenient to express 7 in terms of two parameters 
only, the index of anomalous diffusion 6 — (ji — fl)/v and 
the Hausdorff dimension of the infinite percolation clus- 
ter, df = d — (3/v (56|. Here d is the topological (integer) 
dimension of the ambient space. Remembering that the 
considered system of equations, Eqs. ([T]), is 2D, we have 
d = 2. Nevertheless, the above relation, which is called 
the hyperscaling relation, holds in all d > 2 (5?| ■ In the 
physics context discussed in this work the Hausdorff di- 
mension df measures the subset of phase space occupied 
by the random motions. This dimension is to be distin- 
guished from the fractal dimension of the equipotential 
lines, dfi- The implication is that the particles (random 
walkers) do not follow the equipotential lines exactly be- 
cause of the fine structure of the resonance layer. When 
df is smaller than d, the random motions are, by defi- 
nition, not space-filling. Indeed they fill only a fractal 
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subset of phase space, characterized by df — d — /3/v. 
Eliminating the percolation indices in Eq. © one gets 

7 =(2-d + d/)/(2 + 0). (10) 

To take the non- fractal limit in Eq. (fTTJ|) , one sets the 
Hausdorff dimension equal to the ambient dimension, and 
the index of anomalous diffusion equal to zero. That is, 
df = d and 9 = 0. As a consequence, Bohm scaling 
with 7 = 1 is recovered from Eq. (fTOjl . This limit, be- 
cause of the space-filling property, corresponds with the 
assumptions of classical chaotic behavior. Next we focus 
on fractal (not space- filling) case: df < d. Setting d = 2 
in Eq. (fTO)) one obtains 

7 = d//(2 + 0)=d fl /2. (11) 

Here d s = 2d / 1(2 + 9) is the so-called spectral fractal 
dimension [54], [53] , which has come into physics as the 
density-of-state exponent for vibrational excitations of 
fractal networks J58|. It is noted that d s < df (because 
of 9 > 0). The spectral dimension has the sense of effec- 
tive, fractional number of degrees of freedom in fractal 
geometry [59|. We draw attention to the fact that this 
fractional number is not determined by df only, but also 
requires 9. In general, df and 9 are two independent 
parameters, whose tradeoff defines d s [and hence the 7 
value, according to Eq. (fTTj) ]. The situation is somewhat 
simpler for percolation, as we now proceed to show. 



B. Numerical estimate 

In recent years there has been much excitement about 
the Alexander-Orbach (AO) conjecture [Hf that d s = 
4/3 for d > 2 at percolation. The AO conjecture is exact 
in the high dimensions d > 6 where it holds as a mean- 
field result [60] . In dimensions lower than these, the AO 
conjecture proves to be not exact, being, nevertheless, 
a remarkably accurate estimate of the d s value. Based 
on the available numerical and theoretical predictions, 
the true value of the spectral dimension is believed to 
be slightly smaller than 4/3. The actual discrepancies 
lie within half a percent. See Refs. [54], [13, [Ml- As an 
example, consider the analytic prediction d s = 1.327 ± 
0.001 for d < 6, suggested in Ref. [62|. If, in Eq. (fTT|. 
we apply the mean- field estimate d s =4/3 thus ignoring 
the above small discrepancy we get 7 = 2/3. This value 
coincides with our previous finding (Eq. (28) in Ref. [l2l ]. 
where the Hurst exponent H = 1/2 for diffusion). If one 
is a purist and wants to account for the departure from 
mean- field properties, the procedure is to expedite the 
spectral dimension d s ~ 1.327 [13] or the like 63] to find 
7 ~ 0.66, at little consequence for the rest of the analysis. 
The basic phenomenon is contained in Eq. ©. 



C. Summary 

To this end, our final result for the percolation scaling 
is: D/oj oc Q 2 / 3 , where Q ^> K ^ 1, and we have omit- 
ted the logarithmic correction factor for simplicity. When 
the random motions fill the ambient space (df = d s = d), 
Bohm scaling with 7 = 1 is recovered. A deviation from 
Bohm scaling associated with a fractional 7 smaller than 
1 can be thought of as a signature of concentration of 
the turbulent transport on a subset of phase space with 
fractal geometry. 

IV. FRACTIONAL KINETIC EQUATIONS 
A. Fractional diffusion equation 

The power-law dependence D/uj cx Q 1 which we asso- 
ciate with the initial diffusion on a fractal cluster at per- 
colation must have implications for the generalized form 
of the diffusion equation at time scales a/uj_ <C t < r c <C 
1/lu, for which the charge carriers move only on the frac- 
tal. Observe that the dispersion in Eq. ^ is sub-linear 
at these time scales because of the inverse power-law be- 
havior /(t/r*) oc (i/T*)~' t/,y(2+e) and that the function / 
has not yet reached its asymptotic value. The sub-linear 
dispersion is explained by the trapping effect caused by 
cycles and dead-ends of the fractal (i.e., the "nodes- links- 
blobs" model) (Ref. [57j and Fig. 1 therein). In general, 
obstacles and traps act in a way as to introduce memory 
into the motion. One possible way to include the mem- 
ory is to generalize Fick's first law to (where, to simplify 
the notations, we denote the position vector simply by r) 

j(t,r)=- f A(t-t')Wp(t',v)dt' (12) 
Jo 

with the continuity condition written as 

^p(t, r) = f* A(t - t')W 2 p(t', r)dt'. (13) 

In the above, j(i, r) is the probability current, which flows 
against the concentration gradient, and A(t— t') is a mem- 
ory function, which is nonzero for time lags, for which the 
dynamics are influenced by the fractal geometric proper- 
ties, and is identically zero or vanishing otherwise. In 
order to correctly introduce the memory function, we de- 
fine first a generalized / function via 

i, \ f if s < ,, A . 

M = {/(.) if.>o- (14) 

Thus, by its definition, f(s) coincides with the f(s) func- 
tion for s > and is identically zero for s < 0. Without 
loss in generality, we can set the asymptotic value of f(s) 
to 1, thus exposing the definition of f(s) as follows: 

r if s < 

f(s) = \ ^-yr( 7 ) if < s « s c , (15) 
I 1 if s > s c 
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where s c cx r c is a crossover time scale, where the be- 
havior changes from a power-law to the asymptotic con- 
stant value, and the gamma-function T(-f) is introduced 
for the normalization reasons. Mathematically, the f(s) 
function offers a suitable generalization of the Heaviside 
step function by incorporating an intermediate, power- 
law behavior for < s -C s c . The Heaviside function is 
recovered from Eq. (fT5|) in the limit 7 — *■ 1 . Observe the 
following properties of the f(s) function: (i) /(s) is dis- 
continuous in the origin for all 7 < 1; and (ii) the integral 
J f(s)ds converges for s — > +0 when 7 > 0. Likewise the 
Heaviside step function and Dirac's delta function, the 
f(s) function and its derivatives must be considered as 
generalized functions. Now the memory function is de- 
fined as the time derivative of the f(s) function: 

AGO = ^f(s). (16) 

We have, by means of Eq. (fTr])) . connected the memory 
function with the generalized scaling function /(s) in the 
random walk model. Clearly, A(s) — > for s ^S> s c and 
A(s) = for s < 0. The latter condition is also re- 
quired by causality. When 7 — > 1, the A(s) function 
coincides with Dirac's delta function: A(s) = S(s). In 
this limit, the relation between j(t, r) and Vp(i,r) is lo- 
cal in time [see Eq. (|12p]. Accordingly, the dynamics are 
memoryless, corresponding to an ordinary Fickian diffu- 
sion and the frequency-independent asymptotic diffusion 
coefficient. Combining Eqs. (fT3|) and (fT6|) we have 

m p{t > r) = J* Jt f{t ~ y2/,(i '' r)dt ' (17) 

or 

§- t p(t,r) = - jT ±f(t-t')V 2 p(t',r)dt', (18) 

where the time derivative d/dt' acts only on fit — t'). 
Integrating by parts in Eq. (| 18[) with the initial condition 
V 2 p(i, r) = for t = it is found that 

^p(t,T) = jj{t-t')^p{t',r)dt>, (19) 

where use has been made of /(0) = 0. Focusing on the in- 
termediate, self-similar range, we take f(s) — s 7 ~ 1 /r(7), 
making it possible to rewrite Eq. (fT9|) as 

at^ r) = W)L (^0^aF V p(f ' r) - (20) 

The operator on the right hand side of Eq. (|20|) acting on 
V 2 p(i', r) is known as the Caputo fractional derivative of 
order < 1 - 7 < 1, ^D]' 1 (Ref. for details). With 
this last definition, Eq. (|20|) takes the form of a fractional 
diffusion equation 

?-p{t,r)=lD l t -^p{t,r). (21) 



For 7 — > 1, the Caputo derivative acts as a unit oper- 
ator, thus yielding the familiar — "integer" — diffusion 
equation with no integro-differentiation added. The same 
integer equation is obtained by identifying /(s) with the 
Heaviside step function in Eq. P^|) . 

The end result of the discussion above is that the tur- 
bulent transport associated with a simple random walk 
process on a fractal cluster at percolation is described by 
a fractional time diffusion equation with the fractional 
derivative in the Caputo sense. This equation is differ- 
ent from the equation considered in Refs. [13, 01 
in that it uses the Caputo fractional derivative, o^t" 7 
(instead of more familiar the so-called Riemann-Liouville 
fractional derivative [2l|)- 

B. Fractional relaxation equation 

The characteristic function of the fractional diffusion 
equation, Eq. (|2T|) . obeys the fractional relaxation equa- 
tion with the Caputo time derivative 

^p(*,k) = -k 2 §U t 1 -V(i J k). (22) 

Initial-time behavior of pit, k) can be obtained by reduc- 
ing the fractional derivative to 

^V^k^^^-Vfrk). (23) 

Substituting into Eq. (|2"2"]) and carrying out the trivial 
integration over t leads to a stretched exponential form 
for the p(£,k) function (see, also, Ref. [22j and Eq. (B.7) 
therein): 

p(i,k)~exp{-k 2 tVr(7+l)}- (24) 

This stretched exponential relaxation describes the decay 
of charge density inhomogeneities in self-similar geome- 
try. To obtain a 2D specific expression we set 7 = 2/3 in 
Eq. ([211), yielding p(t,k) ~ exp{-k 2 i 2 / 3 /r(5/3)}. Apart 
from the plasma physics application discussed above we 
expect this relaxation pattern to also characterize charge 
relaxation in thin films of disordered solid materials. 
Some discussion of these properties can be found in Refs. 

Blsilsl. 

The main conclusion to be drawn from the above anal- 
ysis is that the anomalous diffusion associated with the 
random walks at percolation has non-Markovian char- 
acter at intermediate time scales that are in the range 
a/u± <C t < r c <C 1/oj. For t ^> r c , the fractal cor- 
relations, included in qD]~ j , vanish. Consistent with 
this loss of correlation the fractional diffusion and re- 
laxation equations discussed above cross over to their 
integer-derivative counterparts as time t — * 00. It is this 
loss of correlation that permits one to speak about the 
asymptotic transport process in terms of "diffusion," in 
the true sense of the wording, and to obtain the "anoma- 
lous" scaling law in Eq. ([7]) by limiting the power-law 
range of the /(s) variation. 
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V. INCLUDING INERTIA EFFECTS 

The next contribution to the theory of turbulent dif- 
fusion involves the effects of particle inertia in the slowly 
varying electric and magnetic fields. Here, we consider a 
simplified model, which assumes the magnetic field to be 
uniform in space and time: B = const. This model cap- 
tures the essential physics due to inertia. Further gener- 
alization allowing the magnetic field to also be a function 
of time is basically obvious. In a constant B field with a 
slowly varying perpendicular electric field the cross-field 
guiding center drift velocity is accurately approximated 
by (top sign for electrons) [3J| 



uj_(i) = 



EixB 

B 2 



1 dE ± 

u r B dt 



(25) 



where the last term accounts for the polarization drift 
and the time variation of electric field is assumed to be 
slow compared to the gyro frequency. Assuming a peri- 
odic time dependence Ej_ oc e iut with the characteristic 
wave frequency u> <C uj c one finds the velocity compo- 
nents to be 



of the corresponding velocity vector defined as u x = 
u x =F (iuj/u) c )u y and u y = u y ± (iu/u) c )u x , and the char- 



acteristic flow speed u± 

With these definitions one introduces the effective Kubo 
number by means of Q ~ u±/uj^j_. That is, 



(31) 



Accordingly, the effective diffusion coefficient is obtained 
from Eq. ([5]), in which one uses Q instead of Q. Remem- 
bering that 7 = 2/3 at percolation, and keeping the first 
non-vanishing correction due to inertia, we find 

D ~ (In A') 1/3 (a/2a) 2/ Vi(l - J 2 /3^ C 2 )Q 2 / 3 . (32) 

Prospective applications of the results obtained can be 
proposed for the impurity transport in fusion plasma. 
Recently, the turbulent transport of inertial impurities 
have been studied by numerical simulations in Ref. [67j 
and the specific physics consequences arising from com- 
pressibility have been placed in the center of attention. 



E y =F {iu)/u> c )E x , 



(26) 



VI. OVERALL SUMMARY AND FINAL 
REMARKS 



u y = —E x =F {iv/w c )Ey, (27) 

where B = 1 for simplicity. Eliminating E x in Eq. (j26|) 
yields 

u x =F {iLu/uj c )uy = E y - (uj/uj c ) 2 Ey sa E y . (28) 

Likewise in Eq. (|27| : 

u y ± (iu/uj c )u x = -E x + (uj/lo c ) 2 E x w -E x . (29) 

Last terms on the right hand side of Eqs. (|28|) and (|29|) 
could be neglected since they involve the same compo- 
nents of electric field as the leading terms and are smaller 
by (uj/uJc) 2 . Applying d/dx to Eq. (J2SJ); then d/dy to 
Eq. (j29]k and adding the two equations together one finds 

du x du y _ _^_iu^ / du y _ du x \ 
dx dy lu c \ dx dy J 

Here, use has been made of Ej_ = — V$. Equation ([ST)]) 
shows that the polarization drift introduces subtle com- 
pressibility into the motion. Clearly, the compressibility 
effect scales linearly with the mass-to-charge ratio (via 
the dependence on the cyclotron frequency). Topologi- 
cally, the compressibility of advecting flow gives rise to 
limit cycles and stable foci, which attract and trap tracer 
particles (Ref. [ij and Fig. 30 therein). Thus, we expect 
the turbulent diffusivitics of inertial tracers to be compar- 
atively smaller than the diffusivities of ideal particles due 
to compression by the polarization drift. Equations (|28p 
and (f2"9")l suggest the compression effect be described in 
terms of the effective velocity flow, with the components 



In the present work, we have exposed a few crucial 
physics issues behind the so-called "percolation" trans- 
port in low-frequency, electrostatic (anisotropic mag- 
netic) turbulence, basing our investigations on the for- 
malism of Hamiltonian dynamics, random walk models, 
and fractional derivative equations. The central prob- 
lem being addressed is the scaling of the turbulent diffu- 
sion coefficient with the fluctuation strength in the limit 
of slow fluctuation frequencies (large Kubo numbers). 
In this limit, the transport is found to exhibit pseu- 
dochaotic, rather than simply chaotic, properties. "Pseu- 
dochaotic" means random non-chaotic dynamics with 
zero Lyapunov exponents [l?], [H, [Hj]. In our descrip- 
tion, pseudochaos occurs as a consequence of the con- 
centration of the resonant motions on a subset of phase 
space with fractal geometry. Because of the strongly re- 
duced phase space, the dynamics are characterized by 
anomalously slow mixing properties associated with the 
vanishing Kolmogorov-Sinai entropy. 

As a simplified, microscopic model we considered a ran- 
dom walk model on a fractal cluster at percolation ( "the 
unbiased ant in the fractal labyrinth") [53[. Based on 
this model, we found the low-frequency, percolation scal- 
ing of the turbulent diffusion coefficient to be given by 
D/oj oc Q 1 with 7 = 2/3 (here Q 3> 1 is the Kubo num- 
ber), in agreement with the work in Ref. [l2T ] and at 
variance with the prediction 7 = 7/10 in Refs. P, [3, 0). 
When the non-fractal limit is taken in the model (i.e., the 
random motions are thought of as space-filling) , Bohm 
scaling with 7 = 1 is reproduced. In this spirit, a de- 
viation from Bohm scaling, associated with the 7 expo- 
nent smaller than 1, can be interpreted as a signature 
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of concentration of the transport processes on a fractal 
geometry. 

Focusing on the non-Markovian properties of the trans- 
port, by introducing a generalized memory function for 
the random walks at percolation, we have derived a frac- 
tional diffusion equation with the time derivative in the 
Caputo sense (as opposed to a more conventional defi- 
nition of the fractional differentiation in the Riemann- 
Liouville sense). In our study the non-Markovian prop- 
erty has occurred as a consequence of fractality and, 
therefore, has had a purely geometric origin. 

Finally, we have discussed a simple generalization of 
the model treatment described above by taking into ac- 
count finite particle inertia. Inertia enters the model 
equations in the form of the polarization drift and leads 
to a decreased particle diffusivity because of compression 
effects. These results may find further application in de- 
scribing the impurity transport in fusion plasma. 

We have, in the present work, significantly simplified 
the presentations by assuming a characteristic micro- 
scopic time and spatial scales of the random motion. Ex- 
tensions to include a distribution of time intervals be- 
tween consecutive steps are straightforward and have, 
in the case of homogeneous support, discussed in Refs. 
[22I I23I ] on the basis of continuous time random walks 
(CTRWs) [55]. In general we expect the effect of time 
scale distribution to slow down the anomalous diffusion 
due to fractal geometry. In this respect, a fractional 7 
smaller than 2/3 may be conceivable. We remark, how- 
ever, that, while the detailed, microscopic picture of the 
random motion may vary, we expect the basic physical 
properties discussed here to remain essentially the same. 

CTRW-like models and their derivatives are further 
gen eralized to include a distribution of jump lengths 
|22l [H, [2i| [fjl, [69[ physically corresponding to nonlo- 
cal transport in phase space. On the level of CTRWs, 
fractional diffusion models of perturbative transport in 
magnetically confined fusion plasma including nonlocal 
transport have been obtained in Refs. [70l . ITU I72I [z3, 
where one also finds a discussion of the numerical simu- 
lation results. 



There can be various physical mechanisms at play to 
give rise to nonlocal transport. One such mechanism can 
be associated with mode coupling and build-up of cor- 
relations [not included in the wave-like Hamiltonian in 
Eq. ([2])] leading to the formation of large-scale coherent 
structures in the turbulent flow. It is found in direct 
numerical simulations of forced and dissipative turbu- 
lence that the presence of coherent structures leads to 
a spatially nonuniform transport [7J|. Consistent with 
this property, in electrostatic drift-wave turbulence, co- 
herent, vortex structures are found to enhance the diffu- 
sion in the direction of the background density gradient 
[3ll [75j . The propagation of coherent structures can also 
cause significant broadening of the turbulent region and 
in this fashion affect the scaling properties of the trans- 
port [76l Fzrj ] . In burning plasmas, where the energetic 
ions (MeV energies) and charged fusion products con- 
stitute a significant fraction of the total plasma energy 
density, the coupling between the nonlinear energetic 
particle modes [78]], mediated by the energetic particles 
themselves, results in the transition to strong, convective- 
like transport by radially amplifying avalanches [7^, [8(| • 
The process obeys complex nonlinear parabolic equation, 
which, under some nonrestrictive assumptions regarding 
the shape of the energetic particle source function, can 
be cast [28[ in the generic form of a fractional nonlinear 
Schrodinger equation (2?| ■ Indeed the latter equation, 
which is closely related with the fractional Ginzburg- 
Landau equation [8l| . describes the fractional dynamics 
of coupled nonlinear oscillators with long-range interac- 
tion [82j ], After all, we address nonlocal edge phenomena 
in magnetically confined fusion plasma as for instance the 
problem of anomalously fast response in the plasma core 
to a cold pulse edge perturbation [jjl [83[ . Analyses of 
these general phenomena remains to be carried out. 
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